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ABSTRACT 

The velocity dispersion of galaxies on small scales (r ~ lh~ x Mpc), 
cr i2( r ) ; can be estimated from the anisotropy of the galaxy-galaxy correlation 
function in redshift space (Davis & Peebles 1983). We apply this technique 
to "mock-catalogs" extracted from N-body simulations of several different 
variants of Cold Dark Matter dominated cosmological models, including models 
with Cold plus Hot Dark Matter, to obtain results which may be consistently 
compared to similar results from observations. We find a large variation in 
the value of ai2(lh~ l Mpc) in different regions of the same simulation. We 
investigate the effects of removing clusters from the simulations using an 
automated cluster-removing routine, and find that this reduces the sky-variance 
but also reduces the discrimination between models. However, studying o~ 12 
as clusters with different internal velocity dispersions are removed leads to 
interesting information about the amount of power on cluster and subcluster 
scales. We compute the pairwise velocity dispersion directly in order to check 
the Davis-Peebles method, and find agreement of better than 20% in all the 
models studied. We also calculate the mean streaming velocity and the pairwise 
peculiar velocity distribution in the simulations and compare with the models 
used in the Davis-Peebles method. We find that the model for the mean 
streaming velocity may be a substantial source of error in the calculation of 012. 

Subject headings: cosmology - large scale structure of the universe, cosmology - 
dark matter, galaxies - clustering 



1. Introduction 



The velocity dispersion of galaxies on small scales (r ~ lh 1 Mpc), combined with 
cluster abundance data on intermediate scales and the COBE normalization and galaxy 
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peculiar velocity information on large scales, provides a strong constraint on cosmological 
models by constraining the shape of the matter power spectrum. In this paper we 
investigate a method introduced by Peebles (1976; 1980) and Davis & Peebles (1983, 
hereafter DP83), which uses the anisotropy of the redshift-space correlation function to 
determine the pairwise velocity dispersion on small scales. We shall refer to this method as 
the Davis-Peebles method. 

The galaxy-galaxy correlation function £(r) is one of the canonical statistics used 
in studying large scale structure. A related statistic is the redshift-space correlation 
function, £(r p , 7r), which is a function of the components of the separation in redshift 
space perpendicular (r p ) and parallel (7r) to the line of sight. If the correlation function is 
isotropic in real space, it will be anisotropic in redshift space due to the peculiar velocities 
of the galaxies. Hence the degree of anisotropy of £(r p ,7r) is a measure of the moments of 
the peculiar velocity distribution. 

The first moment of the pairwise velocity distribution, v\2(r), is proportional to f2g 6 
if galaxies trace mass and density fluctuations are in the linear regime, where Q is the 
density of matter in units of the critical density at the present epoch. The second moment, 
<7i2, is the velocity dispersion and measures the kinetic energy of the galaxy distribution. 
This quantity has been used in combination with the Cosmic Virial Theorem to estimate 
Q (DP83), although recently some authors have presented arguments that this calculation 
is not only plagued by extremely large uncertainties from cosmic variance (Fisher et al. 
1994b), but that some fundamental assumptions in the usual formulation of the Cosmic 
Virial Theorem may also be incorrect (Bartlett & Blanchard 1995). 

The first calculation of <Ji 2 on a relatively large survey was done by Davis and 
Peebles (DP83). They calculated <7i 2 for the CfAl redshift survey, a survey containing 
1840 redshifts covering 1.83 steradians in the North galactic hemisphere (Huchra et al. 
1983). In another paper (Somerville, Davis, & Primack 1996, hereafter SDP), we present 
a reanalysis of their work, in which we show that <j 12 for this sample is extremely sensitive 
to the way in which corrections for infall into the Virgo cluster are applied, and that 
the value of cr 12 for this survey is dominated by the clusters. The same calculation was 
done on the Southern Sky Redshift Survey(SSRSl) (da Costa et al. 1991), with results 
of a 12 {lh~ l Mpc) = a 12 (l) ~ 300 km/s (SDP; Davis 1988). Recently a 12 has also been 
calculated for the IRAS survey (Fisher et al. 1994b), the CfA2/SSRS2 survey (Marzke 
et al. 1995), and the Perseus-Pisces survey (Guzzo et al. 1995). These calculations have 
shown a range of values of o"i 2 (l) from 272 km/s for SSRS2 to 769 km/s for Perseus- Pisces 
(see Table 1 of SDP for a summary). 

The values of cr 12 usually quoted for simulations are calculated by measuring the 
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dispersion of the pairwise peculiar velocity field directly using the full three-dimension 
position and velocity information for the halos (Davis et al. 1985; Gelb & Bertschinger 
1994; Klypin et al. 1993). In real redshift surveys, not only are there errors introduced by 
edge effects and selection effects, but oyi is extracted by fitting a model to the correlation 
function in redshift space, £(r p ,7r). This quantity is quite noisy especially for samples with 
small numbers of galaxies. In addition the procedure involves a number of assumptions. It is 
a reasonable question, therefore, whether the values from simulations may be meaningfully 
compared to the observational values. Zurek et al. (1994) applied the Davis-Peebles 
method to mock redshift surveys extracted from simulations of a standard Cold Dark 
Matter (CDM) model, and found that this yields a rather large range of values for (7i 2 . So 
far there has not been a comparison of a 12 calculated using the Davis-Peebles method on 
"observed" simulations of different cosmological models. In this paper, we investigate the 
robustness of oyi using mock redshift surveys extracted from several different cold dark 
matter dominated models, and the ability of this statistic to discriminate between such 
models. We estimate the sky-variance and cosmic variance of the statistic and identify 
sources of error in the Davis-Peebles method. We also investigate the effects of removing 
clusters from the samples. 



2. Cosmological Models and Simulations 

The models and simulations are described in detail in Klypin, Nolthenius & Primack 
(1995, henceforth KNP). All of the models studied here have Gaussian initial fluctuations 
with a Harrison-Zel'dovich scale invariant spectrum, density parameter Q = 1 and Hubble 
parameter h = 0.5 (H = lOO/i kms^Mpc- 1 ). 

The models represented are two variants of Cold Dark Matter (CDM) and one variant 
of Cold plus Hot Dark Matter (CHDM). The "standard" Cold Dark Matter model assumes 
that the fraction of the mass in the universe made of baryons is Q b = 0.1 with the rest of 
the mass made up of a species of non-relativistic, dissipationless particle (the cold dark 
matter). There are different ways of normalizing the spectrum — "unbiased" (b=l), which 
assumes that the galaxy fluctuations trace the dark matter fluctuations, and "biased", 
which assumes that perhaps only especially large amplitude fluctuations lead to galaxy 
formation, i.e. ( — ) = b (—) . We have analyzed an "unbiased" (b=l) simulation 

V P ) galaxy VP/ mass 

(which we call CDM1) and a b=1.5 simulation (CDM1.5). In linear theory, this corresponds 
to rms fluctuations of mass in a sphere of 8h _1 Mpc radius as = 1 for the unbiased model, 
(j 8 = 0.667 for the biased model, or a quadrupole in the angular fluctuations of the cosmic 
microwave background of Q ps -norm = 12.8 fiK or Q ps -norm = 8-5/iK respectively (KNP). 
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The value of Q ps - n orm f° r the unbiased model was consistent with the lowest normalization 
quoted in the first-year COBE results (Smoot et al. 1992). It is too low in view of the 
more recent COBE results, which give Q ps -norm — 18/ifT (Gorski et al. 1996), but this 
model is still of considerable interest. The biased model (CDM1.5) is not compatible with 
the COBE measurement of Q ps - n orm', however, it does predict approximately the correct 
number density of clusters (Frenk et al. 1990), and is useful for comparison. 

The Cold plus Hot Dark Matter (CHDM) model assumes that the non-baryonic matter 
is made of a cold component, as before, and a hot component, generally thought to be a 
massive neutrino. In "standard" CHDM, the ratio of cold to hot to baryonic matter is 
fl c /fl u /fl b = 0.6/0.3/0.1 , corresponding to a single massive neutrino with a mass of about 7 
eV. The power spectrum in our model is normalized to a 8 = 0.667, corresponding to b = 1.5 
or Qps-norm = 17 /iK . More recent data and analyses favoring a smaller fraction of the mass 
in neutrinos, VL U ~ 0.2, are summarized in Primack et al. (1995) and Liddle et al. (1996). 
Our preliminary analysis of new simulations with Vt v = 0.2 indicate that oyi in this model 
does not differ substantially from the Q u = 0.3 model on the scales of interest in this paper. 

The simulations were done using a standard Particle-Mesh (PM) code with a 512 3 
force mesh. Each simulation has 256 3 cold particles, and the CHDM simulations have 
two additional sets of 256 3 hot particles with random thermal velocities corresponding 
to Fermi-Dirac distribution (KNP). The size of the computational box is 100 Mpc (i.e. 
50/i _1 Mpc for h = 0.5) and the smallest resolved co-moving scale is ~ lOO/i -1 kpc. 

The initial fluctuations were generated with the same random numbers for CHDM1 
and both CDM simulations. After running the CHDM1 simulation, it was discovered that 
there were two mistakes in the initial conditions (see KNP), which fortunately are in phase 
and largely cancel. A new simulation with the same cosmological parameters but different 
initial conditions (and mistakes corrected) was run — we call this CHDM2. Extensive 
comparison of the two simulations has shown that the power and velocity differences on 
small scales are no more than 5% at z < 7. In addition, there was a statistical fluke (of 
about 10% probability) in the initial conditions for both CDM simulations and CHDM1: 
the amplitude of the longest waves was a factor of 1.3-1.4 larger than that expected for the 
ensemble average, so the power on large scales is approximately a factor of two larger than 
typical. However, this could be considered a compensation for the finite size of the box, 
or in the case of CDM1, as compensation for the low normalization compared to COBE. 
Also, comparing CfAl to the much larger APM survey shows that the CfAl region has 
unusually high power at these scales (Nolthenius, Klypin & Primack 1996 (NKP), Baugh 
& Efstathiou 1993). In any case comparison of CHDM1 and CHDM2, which has a typical 
amount of power, gives us a measure of cosmic variance. 
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In order to compare our simulations with observations, we need to identify objects 
which correspond to galaxies. There is no satisfactory proscription to do this in simulations 
containing only dissipationless dark matter (Summers, Davis, & Evrard 1995). One 
approach is to assume that galaxies form in regions where the dark matter has collapsed 
to a sufficiently high density. However, a well-known problem with this procedure is that 
when the dark matter halos merge, they quickly lose all discernable substructure. This 
results in a halo mass distribution which includes a few very large mass (~ 10 13 — 1O 15 M ) 
halos in the final epoch of the simulations which should probably be associated with the 
dark matter in cluster cores rather than with the halos surrounding single galaxies (Katz 
& White 1993; Gelb & Bertschinger 1994). In simulations with hydrodynamics and gas, 
the baryons can lose energy and form smaller clumps which remain distinct within the 
overmerged dark matter halos (c.f. Evrard, Summers, & Davis 1994; Hernquist & Katz 
1989; Cen & Ostriker 1992). However, associating these baryon clumps with galaxies is still 
an oversimplification as the details of galaxy formation depend on many interdependent 
processes including dissipation, gas dynamics, star formation, and energy feedback from 
supernova. Exactly how these processes affect galaxy formation is not well understood and 
it is not possible to include these effects in N-body simulations on cosmological scales with 
present computing capabilities. Clearly, for the present study of large scale structure, it is 
necessary to attempt to make use of information from dissipationless simulations as best we 
can. Therefore, we have adopted a scheme, informed by results from simulations with gas 
and hydrodynamics (Evrard et al. 1994), for "breaking- up" the overmerged dark matter 
halos. Although this procedure is ad hoc and has many limitations, it does go a step further 
towards enabling a realistic comparison between real redshift surveys and simulations than 
many previous analyses of simulations. 

The procedure for breaking up the halos, assigning luminosities, and forming magnitude 
limited "sky- catalogs" which mimic the CfAl redshift survey is described in detail in 
NKP. First, halos were identified as mesh cells with a sufficiently high dark particle mass 
overdensity Sp/p at the end of the simulation (z — 0). Halos with a mass above a certain 
cutoff were broken up and the fragments were assigned masses according to a Schechter 
distribution. Fragment velocities were chosen randomly from a Gaussian distribution using 
the rms velocity of the nearest neighbors as the dispersion. Luminosities were chosen 
randomly from a Schechter distribution with the same parameters as the CfAl catalog and 
assigned to the halo fragments assuming that higher luminosity corresponded to higher 
1-cell mass. This results in a distribution with the same selection function as the CfAl 
catalog. The broken-up dark matter halos are hereafter somewhat metaphorically referred 
to as "galaxies" . 

To construct the mock CfAl catalogs, six different "home galaxies" were selected for 
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each simulation, in such a way as to mimic the conditions of the CfAl redshift survey: the 
home galaxies were required to lie in an area with local galaxy density within a factor of 
1.5 of the local CfAl density, and to have a Virgo-sized cluster about 20 Mpc away. An 
"observer" was placed on each home galaxy and radial velocities along the line of sight 
were calculated for each galaxy. Catalogs were created with the same angular boundaries 
as the CfAl North survey {b > 40°, 5 > 0°, 1.83 sr), and also all-sky catalogs containing all 
galaxies with \b\ > 10° (10.4 sr). 

In order to check the effectiveness of our break-up procedure, we have verified that the 
galaxy-galaxy correlation function in real space follows a power law all the way down to 
the resolution limit. Without break-up, the correlation function falls below the power law 
on small scales. We have also checked that hierarchical scaling holds for the galaxies; i.e. 
that the volume-averaged three-point function is proportional to the averaged two-point 
function squared. We have found that the reduced skewness S3 = £ 3 (r)/£ 2 (r) 2 i s constant 
to a good approximation from r = 0.5 Mpc to r = 20 Mpc. These results suggest that our 
break-up procedure leads to a halo distribution with the expected clustering properties, at 
least to this order on these scales. 



In this section we briefly describe the method used to extract the pairwise velocity 
dispersion <t 12 from the redshift-space correlation function £(r p , 7r). Readers should refer to 
DP83 and Fisher et al. (1994a, b) for more details of the method. 

The correlation function in redshift space, £(r p , tt), is estimated by counting the number 
of pairs in a bin in r p (separation perpendicular to the line of sight) and n (separation 
parallel to the line of sight). It is normalized by constructing a catalog of Poisson distributed 
points with the same selection function and angular limits as the data, and counting pairs 
between the data and the Poisson catalog: 



where DD is the number of pairs between data and data, and DR is number of pairs 
between the data catalog and a Poisson catalog. The quantities ur and n D are the minimum 
variance weighted densities (see Davis & Huchra 1982) of the Poisson and data catalogs, 
respectively. In practice, we use a large ensemble of Poisson catalogs in order to reduce shot 
noise and to ensure that no bin has zero pair count. 

Let F(w I r) be the distribution function of velocity differences w for pairs of galaxies 
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Method 
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with vector separation r, and /(u>3 | r) the velocity distribution function averaged 
over the directions perpendicular to the line of sight. The first moment of F(w | r) is 
Uv2{r) = — ([v(x) — v(x + r)] ■ f) where the average is number (not volume) weighted. 
This quantity is also known as the mean streaming velocity, and assuming isotropy it is 
a function only of the magnitude of r. Correspondingly, the first moment of f{w^, | r) is 
(w^) = y Vv2(r)/r where y is the component of r along the line of sight. 



A reasonable form for the distribution function f{w^ \ \ 
moments, is 

w 3 (r) - (w 3 (r)) 



parameterized by its 



/(u>3 | r) oc exp 



(2) 



It has been found emperically from studying observations and N-body simulations (Peebles 
1976; Fisher et al. 1994b; Marzke et al. 1995) that on small scales an exponential form 
(n=l) fits the data better than a Gaussian (n=2) or any higher power of the argument. We 
have also tested this assumption for our N-body simulations and find excellent agreement 
with an exponential in all the models up to scales of r ~ 5h~ 1 Mpc. Recently, Sheth 
(1996) gave a derivation of the exponential form for the distribution function using the 



Press-Schechter approach. Adopting this form for f(u>3 
W3 = Ti — y we have: 



and using r 2 



+ y 2 and 



1 + 60>, 7T) = -J= 



dy 



on (r) 



[1 + £(r)] exp 



-y/2 



7T 



Hoy 



H r 



(3) 



An approximation based on self-similar solutions of the BBGKY hierarchy leads to a form 
for vu(r) (Davis & Peebles 1977): 

FH r 



vi2(r 



[1 + (r/r o y 



(4) 



F is an adjustable parameter of the model. The assumption of stable clustering, i.e. that 
the collapse of the cluster is exactly balanced by the Hubble flow, corresponds to F = 1. 
The results presented here are for F = 1 unless stated otherwise. Later on in this paper 
we investigate this model, including the validity of the stable clustering assumption, by 
calculating Uu(r) in the N-body simulations. The velocity dispersion au is then obtained 
by fitting the model of equation (|3p to £(r p , ti) estimated from equation (JTJ) . 



4. Results 



One approach to comparing observations to simulations is to try to translate the space 
of observed quantities (i.e. redshift) into the space of the simulations (real space). Reflecting 
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this approach, DP83 attempted to correct the measured redshifts by modeling the flow 
field around Virgo, which is in the foreground of the CfAl survey. However, we have shown 
that the value of ayi depends sensitively on the details of these corrections (SDP). Another 
approach is to "observe" the simulations to form simulated redshift surveys, which attempt 
to mimic as closely as possible what astronomers would actually observe if they lived in 
the universe of that simulation, and to then analyze the simulated and observed catalogs in 
exactly the same way Taking this approach, we have selected catalogs from the simulations 
requiring that a Virgo-sized cluster appear in the foreground as described in section |, and 
do not include corrections for cluster infall in any part of the procedure. 

Contour plots of £(V p ,7r) for the 10.4 sr sky-catalogs are shown in Figure For 
each simulation, £(r p ,7r) has been averaged over 6 views and smoothed. (For purposes 
of the contour plots only, £(r p , n) was calculated with linear bins in r p and n. For the 
rest of the calculation, logarithmic bins were used in r p and linear bins in ir). There is a 
visible difference in the shape of the contours for the different models. In the absence of 
peculiar velocities, the contours would be perfectly circular segments. However, the velocity 
dispersion of clusters causes structures to appear elongated along the line of sight for small 
r p (the familiar "fingers of god"), and infall causes the contours to appear compressed in 
the line of sight direction for larger r p . This effect is much more pronounced in the CDM 
models as expected due to the larger velocity dispersions and infall velocities. 

The values of and 7 obtained from the inversion of w(r p ) (see DP83) and used in the 
fit for a 12 are given in Table 3. The results for <7i2(r) are shown in Figure |2|. The points are 
the average of the results for the six views, and the error bars are standard deviations for 
the different views. These error bars reflect the "sky variance" , or variation within the same 
simulation when viewed from different points, which for this relatively small box size is 
likely to be an underestimate of the true cosmic variance. The difference between CHDM1 
and CHDM2, about 200 km/s at lh~ x Mpc, may give a better estimate of the possible 
cosmic variance. This is consistent with the large variance seen in recent calculations of a±2 
in different redshift surveys (see SDP). 



4.1. Removing Clusters 

We attempted to address this apparent non-robustness of au by developing an 
algorithm to automatically remove clusters from catalogs in a way which could be applied 
consistently to both real and simulated data. It is no surprise that rich clusters are a major 
source of sample-to-sample variation in cr 12 . Because £(r p ,7r) is a pair-weighted statistic, 
clusters tend to dominate it. We saw in SDP that removing the Virgo cluster from the 
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CFAl survey significantly reduced 012, consistent with the findings of Zurek et al. (1994). 
The SSRS1 survey does not contain any clusters as rich as the Virgo or Coma clusters found 
in CfAl, and oyi is much smaller for SSRS1 (0x 2 (l/i _1 Mpc) ~ 320 km/s) than for CfAl 
(o-i 2 (l/i _1 Mpc) ~ 620 km/s)(SDP). Fisher et al. (1994b) found 12 (1) for the IRAS redshift 
survey to be only about 320 km/s, much lower than the result for CfAl. The IRAS survey 
is dominated by dusty spiral galaxies and undercounts cluster centers by about a factor of 2 
relative to optically selected surveys. Marzke et al. (1995) also found that removing Abell 
clusters of richness class R>1 significantly reduces 0x2 in the CFA2/SSRS2 survey, and that 
a 12 changes more drastically in regions of the survey where there were many rich clusters to 
begin with. In view of this evidence, we thought it interesting to see what happens to oyi 
when the clusters are removed from the simulations. We therefore developed a method to 
remove clusters which can be applied consistently to both simulations and observations. 

Our algorithm is as follows. We divide the catalog into bins and identify bins with 
density fluctuations larger than a specified cutoff. We then calculate the luminosity 
weighted centroid of all points lying within a cylinder of radius r c and redshift interval 2h c 
centered on the bin where the fluctuation was found. We use h c = 1000 km/s, and adjust r c 
as described later. We then take a new cylinder centered on the centroid, calculate the new 
centroid, and continue to iterate in this way until the position of the centroid changes by 
less than some small value. Finally, we calculate the velocity dispersion of all the galaxies 
lying within the cylinder around the converged centroid, and cut all the galaxies in this 
cylinder only if the velocity dispersion is greater than a cutoff a c . We have tested our 
algorithm by visualizing the surveys to make sure that the regions which are cut correspond 
to those that would be identified visually. The number of clusters identified depends on 
both parameters r c and a c . We chose a value of r c by plotting profiles of the clusters and 
identifying the radius at which the number density had dropped to the background level. 
The fiducial value we chose, r c = 2.0/i^ 1 Mpc, is fairly close to the usual Abell radius. 
Table 1 shows the results of varying a c . As a c is lowered, the algorithm identifies more 
and more objects as "clusters" — of course as we go to lower values of a c we are really 
starting to identify objects which we would normally refer to as groups. We select o c = 500 
km/s as corresponding to what are usually referred to as clusters because at this cutoff 
approximately 4-9 % of the galaxies are in clusters, which roughly agrees with observations. 
Table 2 shows properties of clusters identified at o c = 500 km/s. 

The results for 0x2(7") after cluster removal with r c = 2.0h~ 1 Mpc and a c = 500 km/s 
are shown in Figure The sky-variance has decreased somewhat for the CDM models 
but not significantly for the CHDM models. The cosmic variance between the two CHDM 
models has decreased slightly: 012(1) for CHDM1 and CHDM2 now differs by about 150 
km/s instead of 200 km/s. However 012(1) for the biased CDM model now lies between 
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the two values for CHDM1 and CHDM2, and even standard CDM gives a value of 012(1) 
which is within the view-to-view error bars of the CHDM1 value. Therefore, the small 
improvement in robustness acheived by removing the clusters appears to have been attained 
at the expense of discrimination between models. This problem cannot be solved by using 
a different value for a c , as can be seen from Figure f|. The discrimination becomes worse as 
we go to lower values of a c , and for higher values there is very little change in the values 
of cr 12 . The fact that 012(1) appears to converge to almost the same value for all of the 
models as we lower a c is probably an indication that, once we remove the collapsed objects, 
012 is really a measure of Qq. In fact, originally this statistic was not designed as a tool 
to discriminate between models. It was hoped that it would be useful as a measure of Qq 
(or of (3 = f(Q)/b, to be precise, where b is the bias factor defined in section Q) (DP83, 
Davis 1995). It would be interesting to see what would happen as we remove clusters in 
simulations with different values of Qq- 

The change in 012 as a function of a c has a different shape for the different models. 
In Figure |5|, we show the change in 012(1) as a function of the number of clusters cut at 
different values of a c . It is interesting that the curve for CHDM1 lies on top of the curve 
for CHDM2 on this plot even though the values of 012 are very different. What this seems 
to indicate is that although the curve for CHDM1 extends further to the right, indicating 
that CHDM1 has more clusters than CHDM2 (because of the excess large scale power), 
the clusters themselves have the same velocity structure, because the small scale power is 
the same. The CDM models have more power on small scales and the clusters have more 
kinetic energy, as indicated by the steeper rise in A0i 2 (l). In addition the curves extend 
further to the right than CHDM2 because the CDM models also have a larger number 
density of clusters than a typical CHDM model. Plotted in this way, this quantity may give 
interesting information on the amount of power present on cluster and subcluster scales. It 
appears to be more robust and discriminatory than the actual value of 012(1)- However, 
it is not very useful at present because existing redshift surveys do not contain enough 
clusters to allow this quantity to be calculated in the region where it is discriminatory. It 
would be interesting to study this in future large volume surveys with good sky coverage. 

4.2. The Form of the Pairwise Peculiar Velocity Distribution and the Mean 

Streaming. 

A form for the pairwise peculiar velocity distribution must be assumed in the procedure 
we have used to calculate 012. As we discussed briefly, the form used by DP83 and other 
authors who have recently completed similar analyses is an isotropic exponential, which 
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was proposed by Peebles (1976) because it appeared to fit observations. This functional 
form was further investigated by Fisher et al. (1994b) and Marzke et al. (1995) and found 
to be consistent with the IRAS and CFA2/SSRS2 redshift surveys, and with Cold Dark 
Matter N-body simulations. However, we thought it would be interesting to study the form 
of the pairwise velocity distribution in real space for our different models to re-evaluate 
whether this is the most appropriate form. We also take advantage of the fact that unlike 
in real observations, in the simulations we know the real-space positions of the galaxies and 
their peculiar velocities separately. We compare the "true" value of au obtained by fitting 
directly to the peculiar velocity distribution with the results of the Davis-Peebles procedure. 
We may then evaluate whether the "true" value of au, uncomplicated by extracting it from 
the redshift space data, is a robust discriminator between models. In this way we can test 
the accuracy of the Davis- Peebles method for extracting <7i2(r). 

We find that the exponential is an excellent fit to the distribution on small scales in all 
the models. On scales of r ~ bh^ 1 Mpc the distribution begins to look a bit flatter at small 
w 3 than the exponential, and the distribution begins to approach a Gaussian at even larger 
separations. However, at the scales where this procedure is used (~ lh^ 1 Mpc), using this 
model is not likely to be a significant source of error. 

As we mentioned in section [| following DP83 (and subsequent workers) we have 
modeled the mean streaming velocity v±2(r) by equation @ in our fitting procedure. 
Marzke et al. (1995) and others have suggested that the observed scale dependence of o\i 
may be merely an artifact of this term. We have investigated this quantity by computing 
it directly in the simulations. It is shown in Figure |] along with the model for F=l and 
F=0.5. The general form of the model holds on intermediate scales, but there is a large 
amount of variation between the different views, and the prediction of the model on scales 
of Xh~ x Mpc is in many cases quite inaccurate. The CDM models show a rather large 
negative mean streaming velocity on small scales, which may correspond to shell crossing. 
On average the model with F = 1 appears to overpredict W12 for CHDM, while F = 0.5 
may be a better fit. This is surprising because F = 1 corresponds to streaming which 
exactly cancels the Hubble flow on scales less than the correlation length (stable clustering), 
F > 1 corresponds to collapse on those scales, and F < 1 indicates that the clusters are not 
collapsing, and in fact are expanding. This implies that stable clustering on is not a good 
assumption for the CHDM models. 

To show the dependence of the results on the assumed value of F, in Figure [7] we 
show the results obtained for au (using the Davis-Peebles method on the "observed" 
sky-catalogs, as before) when we do not include the streaming model in the fit, i.e. when 
we set F = 0. This systematically reduces the values of au that we obtain, and changes 
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the scale dependence as expected. The CHDM models are now entirely flat over r p . It is 
also possible to allow F to be fit as a free parameter. The problem with this approach is 
that there is a degeneracy between F and (r 12 — we tended to obtain lower values of both 
parameters when we allowed F to be fit freely. In addition, having two free parameters 
increases the error on the fit. Apparently, the use of this model for the mean streaming, 
especially with the assumption of stable clustering (F — 1) could be a substantial source of 
error in the procedure. 

Figure || shows the results obtained when i2 was fit directly to f{w^\r) using the 
peculiar velocity information from the simulations. Although we use the real space 
positions and peculiar velocities of the galaxies, we use the same galaxies as those in the 
magnitude-limited sky-catalogs in order to facilitate a direct comparison of <J trU e with the 
value obtained by the Davis-Peebles method (a dp)- Table 3 shows odp with and without 
the streaming model and <Jt rU e at lh" 1 Mpc, averaged over the six views. The results for 
a dp are within 20% of atrue for all the simulations. The agreement with Otme is better for 
all the simulations except CHDM2 if we do not include the streaming model {F — 0). The 
view-to-view variation is comparable for odp and a true , which suggests that the variation is 
not noise introduced by the Davis-Peebles method of extracting i2 but that it is intrinsic 
in the simulations and presumably in the real Universe. 



5. Discussion and Conclusions 

For the N-body simulations of three different cosmological models that we have 
analyzed, we find that the values of a 12 are considerably higher for CDM models than for 
CHDM models. The CHDM models, which give 012(1) ~ 540(440) km/s for CHDM2 to 
740 (600) km/s for CHDM1 are perhaps more consistent with the body of observational 
values taken as a whole. (Numbers in paranthesis are for F — 0). Although 012(1) = 647 
km/s for CfA2 North is marginally consistent even with unbiased CDM (012(1) ~ 1024(880) 
km/s), 0i 2 (1) for CfA2 South (367 km/s) and SSRS2 (272 km/s) are not (Marzke et al. 
1995). However, because of the problems with galaxy identification discussed in Section 
|2], we do not think that any model studied here should be ruled out on the basis of these 
results. Any existing large volume N-body simulations would have similar problems. The 
point of this paper is precisely to suggest that it is premature to use cr 12 to draw any strong 
conclusions about cosmological models. However, studying the simulations has given us 
other interesting information. 

We have estimated the expected sky-variance and cosmic variance of 012 in our models. 
The following values are quoted for 0i 2 (l/i _1 Mpc) . The sky- variance (calculated as the 
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standard deviation over six mock catalogs) ranges from ~ 40 km/s for CHDM2 to ~ 145 
km/s for CDM1. The cosmic variance (between two simulations with different initial 
conditions) for CHDM is ~ 200 km/s. The errors usually quoted for <t 12 are formal errors 
on the fit (for which we obtain typically ~ 40 km/s) and in general are underestimates 
of the actual statistical errors. We evaluate the accuracy of the Davis-Peebles method 
for extracting (j\2 from redshift catalogs by comparing with the results of computing the 
velocity dispersion directly in our mock catalogs. We obtain agreement of better than 20% 
for all of our models. This leads us to interpret the large range of values of a 12 obtained 
from different redshift surveys as an intrinsic variation due to the sensitivity of the statistic 
to the clusters contained in the sample, rather than being due to errors in the method. 

We have investigated the effects of removing clusters from the samples using an 
automated procedure. It was hoped that this might make ayi a more robust statistic. 
However, we found that although this reduces the sample-to-sample variation in ayiij) by 
a small amount, it actually reduces the ability of the statistic to discriminate between the 
cosmological models we studied. This may be due to the fact that all of our simulations are 
of Q = 1 models, and that once clusters are removed ayi is really a measure of However, 
our study of the simulations suggests that the change in 0"i2(l) as a function of the number 
of clusters removed may be an interesting quantity to study in future redshift surveys. 

We find that an exponential form for the pairwise peculiar velocity distribution is 
an excellent approximation on small scales (r < Mpc) in all the models studied. 

Measuring the mean streaming Vu{ r ) directly from the simulations revealed that although 
the general form of the model used in the Davis-Peebles method does hold, on small scales 
the measured values may deviate from it considerably. We found that stable clustering is a 
reasonable approximation in the unbiased (b=l) CDM model but not in the CHDM models. 
The use of the BBGKY model for the mean streaming, especially with the assumption of 
stable clustering (F = 1) could be a substantial source of error in the Davis-Peebles method. 

We have shown that o\2 is very sensitive to both the number and the properties of 
the clusters in a sample. This makes <7i 2 a poor constraint on cosmological models given 
the current situation with regards to both simulations and observations. Even the largest 
existing redshift surveys do not represent a fair sample of rich clusters. Also, current N-body 
simulations do not simulate clusters realistically because cluster properties are probably 
sensitive to non-gravitational physics such as gas hydrodynamics, star formation, and 
supernova feedback — effects which are impossible to include in large volume simulations 
with current computing capabilities. As larger redshift surveys become available and it 
becomes possible to simulate clusters in a cosmologically relevant volume, perhaps the 
robustness of 0"i2(l) will improve. In fact, we have suggested a way in which the very 
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sensitivity of a 12 to the properties of clusters could be used to define an interesting statistic 
for characterizing large scale structure in larger samples. In the meantime it is worthwhile 
to work on developing statistics which are discriminatory but less sensitive to the properties 
of clusters. 
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a cut (km/s) 


N c i 


Ngal 


fcl 


n cl {10- 7 Mp C - 3 ) 


CDM1 


600 


8 


106 


1.2 


1.5 




500 


37 


773 


8.8 


6.9 




400 


78 


1605 


18 


14 




300 


116 


2254 


26 


22 


CDM1.5 


600 


13 


204 


2.3 


2.4 




500 


45 


776 


8.7 


8.3 




400 


81 


1398 


16 


15 




300 


121 


2142 


24 


22 


CHDM1 


600 


5 


89 


1.0 


.93 




500 


25 


530 


6.0 


4.6 




400 


78 


1738 


20 


14 




300 


130 


2865 


33 


24 


CHDM2 


600 
















500 


15 


368 


4.3 


2.8 




400 


34 


800 


9.4 


6.3 




300 


82 


1807 


21 


15 



Table 1: Properties of clusters identified using the algorithm described in the text, averaged 
over six views, for the 10.4 sr sky-catalogs. N c i is the number of clusters with internal velocity 
dispersion greater than a cu t, N ga i is the total number of galaxies in those clusters, f c \ is the 
percentage of galaxies in clusters, and n c i is the number density of clusters. 
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Simulation a c i{km/s) N ga i cri2(l) 

CDM1 550 ± 49 20 ± 16 745 ± 114 

CDM1.5 562 ±45 17 ±9 605 ± 38 

CHDM1 552 ±46 21 ± 10 652 ± 37 

CHDM2 533 ±23 25 ± 22 508 ± 41 

Table 2: Properties of clusters with a cut > 500 km/s. The second column shows a c i, the 
mean internal line-of-sight velocity dispersion of the clusters, with the standard deviation 
shown as the error. Column 3 shows the mean and standard deviation of the number of 
galaxies per cluster. Column 4 is (712(1) after the clusters were removed. 



Simulation cr DP (F = 1) a DP (no streaming, F=0) cr true 

CDM1 1024 ± 145 878 ± 133 862 ± 33 

CDM1.5 764 ± 71 658 ± 62 677 ± 41 

CHDM1 736 ±43 607 ± 40 656 ± 52 

CHDM2 537 ±40 444 ± 35 543 ± 72 



Table 3: Comparison of a u computed using the method of Davis and Peebles with a trU e, 
computed by fitting directly to the pairwise peculiar velocity distribution. 
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Fig. 1. — Contours show £(r p , n) averaged over six views and smoothed. Contour spacings 
are A£ = 0.1 for £ < 1 and Alog£ = 0.1 for £ > 1, with the solid contour indicating £ = 1. 

Fig. 2. — The velocity dispersion for the 10.4 sr sky-catalogs, averaged over six views, with 
the standard deviation shown as error bars. The text markers CfA2N, CfA2S, and SSRS2 
show the values from Marzke et al. (1995), for the CfA2/SSRS2 survey. They are shown 
offset to r p — 1 for visibility but actually were calculated for exactly the same bin in r p as 
our result. 

Fig. 3. — The velocity dispersion for the 10.4 sr sky-catalogs with clusters removed, averaged 
over six views, with standard deviations shown as error bars. The text markers CfA2N, 
CfA2S, and SSRS2 show results from Marzke et al. (1995) for CfA2/SSRS2 with Abell 
clusters with richness R > 1 removed. See comment with Figure Q on the horizontal offset. 

Fig. 4. — The velocity dispersion 012(1), with clusters removed at different internal cluster 
velocity dispersion thresholds (a cut ). These values are for one mock catalog (view a). 

Fig. 5. — The change in velocity dispersion at r = lh^ 1 Mpc (Acr^l) = 012(1) — 
&i2(l-) dusters cut) as clusters are removed at different thresholds, as a function of the number 
of clusters cut (Nd). These results are for one mock catalog (view a). 

Fig. 6. — The first moment of the pairwise peculiar velocity distribution, or mean streaming 
1^2 (r). The dotted line is the BBGKY model (equation ^) with F — 1, the dashed line is 
the model with F=0.5, and the smooth solid line is the Hubble flow (Hq = 50 km/s/Mpc). 
Results are shown for all six views in each simulation. 

Fig. 7. — The velocity dispersion for the simulations with F=0 used in the streaming model 
(equation ^), averaged over six views, with standard deviations shown as error bars. The 
text markers CfA2N, CfA2S, and SSRS2 show results with F=0 from Marzke et al. (1995) 
for CfA2/SSRS. See comment with Figure |^ on the horizontal offset. 

Fig. 8. — The second moment of the pairwise peculiar velocity distribution f(ws | r), 
obtained by fitting the exponential model directly to f(w 3 | r). 
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